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Numerical simulations of line-profile variation 
beyond a single-surface approximation for oscillations in roAp stars 
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Abstract 

Prior to the last decade, most observations of roAp stars have concerned the light variations. Recently some 
new, striking results of spectroscopic observations with high time resolution, high spectral dispersion, and 
a high signal-to-noise ratio became available. Since the oscillations found in roAp stars are high overtones, 
the vertical wavelengths of the oscillations are so short that the amplitude and phase of the variation of each 
spectroscopic line are highly dependent on the level of the line profile. Hence, analyses of the variation 
of the spectroscopic lines of roAp stars potentially provide us with new information about the vertical 
structure of the atmosphere of these stars. In order to extract such information, a numerical simulation of 
the line-profile variation beyond a single-surface approximation is necessary. We carried out a numerical 
simulation of line-profile variation by taking account of the finite thickness of the line forming layer. We 
demonstrated how effective this treatment is, by comparing the simulation with the observed line profiles. 
Key words: line: profiles — stars: chemically peculiar — stars: oscillations — stars: variables: other 



1. Introduction 

A fraction of A-type stars, called Ap stars with the sub- 
script "p" , indicating "peculiar" , show strong absorption 
lines of unusual elements. Spectral analyses of these stars 
indicate that rare earth elements are very abundant, at 
least at a certain atmospheric level. Such anomalous spec- 
tra seen in Ap stars are thought to be a consequence of the 
levitation and resultant accumulation of those elements 
into the high atmosphere due to radiation pressure. The 
photon momentum is selectively transferred to the ions of 
those elements through their appropriate absorption lines. 
Many of the Ap stars have global magnetic fields, of which 
the apparent strength is on the order of 10 mT (= lkG), 
and varies periodically and mostly sinusoidally on a time 
scale of days to decades. This periodic variation is likely 
to be caused by rotation of the star in which the magnetic 
axis is inclined to the rotational axis. The spectrum and 
mean brightness vary synchronously with the magnetic 
field strength. Such synchronous variation is thought to 
be a consequence of the fact that the atmospheres of mag- 
netic Ap stars are laterally inhomogencous, so that some 
chemical elements are selectively levitated near the mag- 
netic polar region, while others are levitated near the mag- 
netic equatorial zone, because the presence of magnetic 
fields strongly influence diffusion. This phenomenological 
model is known as the "oblique rotator model" , and it is 
now widely accepted as explanating the magnetic, spectral 
and the mean light variations seen in Ap stars. 

A-type stars occupy the region of the Hertzsprung- 
Russell (HR) diagram where the Cepheid instability strip 
crosses the zero-age main sequence. The A-type stars lo- 
cated in the instability strip are thus expected to pulsate, 
and there are indeed many pulsating stars - the 5 Scuti 



type variables - whose pulsation period is on the order of 
two hours. Pulsational variability has long been sought 
in Ap stars. It is almost thirty years ago since rapid light 
variations in the magnetic holmium star HD 101065 were 
discovered by Kurtz (1978). Similar rapid light variations 
with periods in the range of 6 to 12 min were discovered 
soon thereafter in several cool magnetic Ap stars (Kurtz 
1982), and the concept of "rapidly oscillating Ap stars" 
was established. Currently, 40 Ap stars are known to pul- 
sate. The period range of the oscillations found in Ap 
stars is 6-21 min and is substantially shorter than the 
dynamical timescale, 2tt[R 3 /(GM)] 1 ^ 2 , which for main- 
sequence A-type stars is typically about two hours. This 
is obviously the reason why the stars are called rapidly os- 
cillating Ap stars (Kurtz 1982), or "roAp stars" for short 
(Matthews et al. 1987). This means that the oscillations 
in Ap stars are of very high overtones. 

Another noticeable characteristic of oscillations in Ap 
stars is that, in many cases, the pulsation amplitudes 
of light variations vary synchronously with the apparent 
magnetic field strength of the star. In some cases, even 
though no magnetic field has yet been detected, the pul- 
sation amplitudes are found to vary periodically with a 
time scale of several days. Kurtz (1982) showed that am- 
plitude modulation of the rapid oscillations in Ap stars is 
well-explained in terms of an axisymmetric dipole mode 
whose axis is aligned with the magnetic axis and inclined 
to the rotational axis of the star. It is then obvious that 
the rapid oscillations in Ap stars are strongly influenced 
by stellar magnetism; in this sense, the roAp stars are 
unique among the various types of pulsating star. Recent 
long time-series photometry has enabled us to perform im- 
proved asteroseismology of this class of stars (Kurtz et al. 
2005b; Cameron et al. 2006; Gruberbauer et al. 2008; 
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Huber et al. 2008; Balona et al. 2011a, 2011b; Kurtz et 
al. 2011). 

As mentioned above, the atmospheres of Ap stars arc 
laterally chemically inhomogencous, offering the potential 
to map the oscillations. Furthermore, since the oscilla- 
tions found in roAp stars are high overtones, the vertical 
wavelengths of the oscillations are so short that the am- 
plitude and phase of the oscillations are highly dependent 
on the atmospheric level. Hence, the analyses of rapid os- 
cillations in roAp stars potentially provide us with a new 
tool to diagnose the inhomogencous — both laterally and 
vertically — atmospheres of Ap stars, which are governed 
by radiative levitation and strong magnetic fields. High- 
resolution spectroscopic observations are required for this 
purpose. 

Prior to the last decade, most observations of roAp 
stars have concerned the light variations [see reviews; e.g., 
Kurtz k Martinez (2000)]. During the early phase of 
spectroscopic investigations, observations were done with 
low-resolution spectroscopy. In the case of a Cir, it was 
found that some lines are apparently pulsating in anti- 
phase with others. Baldry et al. (1998) interpreted this 
in terms of a high-overtone standing wave with a velocity 
node in the atmosphere of the star. Bisector measure- 
ments of the H a line show that the velocity amplitude 
and the phase of the principal oscillation mode vary sig- 
nificantly, depending on height in the Ha line. This fact 
confirms the existence of a radial node in the atmosphere 
of the star. Therefore, it has been suggested that the use 
of spectral lines that form at different levels in the atmo- 
sphere can be used to gain depth information about the 
eigenfunction. 

Thanks to big telescopes and high resolution spec- 
trographs, recently some new, striking results of spec- 
troscopic observations of oscillations of roAp stars with 
high time resolution, high spectral dispersion, and a 
high signal-to-noise ratio became available (Savanov, 
Malanushenko k Ryabchikova 1999; Kochukhov k 
Ryabchikova 2001a, 2001b; Balona & Zima 2002; 
Kochukhov et al. 2002; Balona 2002; Balona k 
Laney 2003; Kurtz, Elkin & Mathys 2003, 2005a, 2005c; 
Mkrtichian, Hatzcs k Kanaan 2003; Sachkov ct al. 2004; 
Elkin, Kurtz k Mathys 2005; Elkin et al. 2005; Hatzes 
& Mkrtichian 2005; Mkrtichian & Hatzes 2005; Kurtz, 
Elkin k Mathys 2006; Kochukhov et al. 2007; Kurtz et 
al. 2007; Ryabchikova et al. 2007a, 2007b; Elkin, Kurtz 
k Mathys 2008; Mkrtichian et al. 2008; Sachkov et al. 
2008; Freyhammer et al. 2008, 2009; Kochukhov ct al. 
2008, 2009). These time-series spectroscopy observations 
clearly show that the amplitude and phase of variation of 
each spectroscopic line are indeed dependent on the in- 
tensity level of the line profile. Hence, the variation of 
the spectroscopic lines of roAp stars is a very promising 
information source about the vertical structure of the at- 
mosphere of these stars as well as the magneto-acoustic 
oscillations in these stars. 

Modeling of line-profile variations (LPV) has so far 
been developed mainly to describe the observed variations 
among O, B variables, — (3 Cephei stars, SPB stars, ( 



Ophiuchi variables, which arc nonradial pulsators oscillat- 
ing with low order p-, or g-modes (Osaki 1971; Stamford 
k Watson 1976, 1977; Smith 1977; Kubiak 1978; Vogt 
k Penrod 1983; Balona 1986a, 1986b, 1987; Kambe & 
Osaki 1988; Aerts k Waelkens 1993; Schrijvers et al. 
1997; Townsend 1997). In these cases, contrary to the 
case of roAp stars, the vertical wavelengths of the oscil- 
lations are much longer than the thickness of the line- 
forming layer of each spectral line, and hence the line- 
forming layer at a given longitude and latitude on the 
stellar surface moves up and down uniformly throughout 
the layer. In order to sec the influence of such oscillations 
on the spectral line profiles, we do not have to model the 
line shape in detail. We may ignore the thickness of the 
line-forming layer, and may treat the layer as a single sur- 
face to take account of the motion in the Doppler shift 
of the spectral line. Actually, in the standard procedure 
of line-profile modeling of these stars, it is sufficient to 
approximate the intrinsic line shape as a Gaussian profile 
to see the Doppler effect of oscillations. Beside that, the 
symmetry axis of their oscillations coincides with the ro- 
tation axis of the star, and hence the aspect angle of the 
oscillation axis does not change at all with time. 

However, the situation is different in the case of roAp 
stars, of which oscillations have vertical wavelengths as 
short as the thickness of the line-forming layers. In order 
to extract information about the structure of the atmo- 
sphere from spectroscopic observations of oscillations with 
high time resolution, a numerical simulation of line-profile 
variation beyond a single surface approximation is neces- 
sary. Furthermore, the aspect angle of the symmetry axis 
of the oscillation patterns changes with time along with 
rotation of the star. Therefore, we need to model prop- 
erly the line profiles of roAp stars to extract asteroseismic 
information from the spectroscopic analyses of the time 
series of line-profile variation of these stars (Shibahashi 
2008). 

The main purpose of the present work is to carry out 
numerical simulations of line-profile variation by taking 
account of the finite thickness of the line-forming layer. 
We demonstrated how effective this treatment is, by com- 
paring the simulation with the observed line profiles. The 
procedure of calculating the intrinsic line profiles is out- 
lined in section 2. The equation of radiative transfer was 
numerically solved. The recipes for computation of the 
absorption coefficients are described in section 3. We will 
verify the numerical code by comparing the computed line 
profile with observations in section 4, and demonstrate the 
difference in LPV caused by motion with different verti- 
cal wavelengths. More detailed calculations are described 
in section 5 based on a realistic model of the A-type star. 
The variation in line bisectors is discussed in section 6, and 
a comparison of the simulation results with observations 
is made. Future prospects are discussed in section 7. 
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2. Procedure of calculation of intrinsic line pro- 
file 

2.1. Outline 

It may be instructive to outline here a procedure to cal- 
culate the intrinsic line profile. We divide the stellar disk 
image by a mesh into small square cells, which are labelled 
by the (x,y)-coordinates. The energy flux going from each 
cell toward the observer per second in a unit wavelength 
at the wavelength A is calculated from specific intensity 
of radiation, I Xl which is described by the equation of 
radiative transfer, 



dl\(x,y,z) 



- («a + h + cf\)ph + jxp, 



(1) 



where z denotes the path length, and p is the density 
in mass per unit volume; k x and t\ mean the continuous 
absorption coefficient per mass and the line absorption co- 
efficient per mass, respectively, <j\ denotes the scattering 
coefficient, and j x is the total emission coefficient. Since 
the beam path of light is different from cell to cell, 7a is 
regarded as being a function of (x, y) and the derivative 
is taken as a partial one in equation (1) (see figure 1). By 
dividing both sides of (1) by — (k\ +l\ + cr x )p, we can 
rewrite the equation of radiative transfer in the form of 

dI\(x,y,T Xz ) 



dr Xz 



I\ — S\ 



where S x is the source function defined by 

S\ = j x /(n x +£ x + a x ) 
and 



(2) 



(3) 



dr Xz = -{n\+£\ + <J\)pdz. (4) 

With (Lt Xz , we introduce the optical depth, t Xz , looking 
from the outside in, along the path for the outgoing beam 
of light. The intensity at the observer's site, I x (x,y,Q), is 
derived by an integral of equation (2) with respect to t\ z 
from to oo, and is expressed as 



h(x,y,0) 



S\{T\ z )exp(-T Xz )dT Xz . 



(5) 



The total energy flux going outward toward the observer 
per second in a unit wavelength band is then obtained by 
an integral of l x (x,y,0) over the stellar disk, 



F x oc l x (x,y,0)dxdy 

./disk 



S\{T\ z )cxv(-T\ z )dxdydT\ z 



(6) 



In order to theoretically determine the profile of a spec- 
tral line for a given model atmosphere, we need to estimate 
the background, continuum, radiative flux from the star at 
the wavelength of the spectral line of interest and to esti- 
mate the flux at the same wavelength by taking account of 
both the line absorption and continuum absorption. The 
former is calculated by setting £ x = 0. Then, the line pro- 
file against the background continuum spectrum is given 
by residual flux, defined as the ratio of these two fluxes. 



c 




Ray path 




observer 



Fig. 1. Stellar disk image subdivided by a mesh coordinated 
with (x,y) into small square cells. At each cell, the energy 
flux going toward the observer per second in a unit wavelength 
is calculated by integrating the specific intensity of radiation 
along the path in the direction of z-axis, which is toward the 
observer. 

2.2. Source function 

The total emission, j x , is composed of the thermal emis- 
sion and the contribution of photons coming there by scat- 
tering from other directions. As for the former, we assume 
here for the sake of simplicity local thermodynamic equi- 
librium (LTE) . Then, the thermal emission is described by 
the Planck function, B x . As for the latter, we assume co- 
herent scattering, which does not change the wavelength 
of photons after scattering, and also isotropic scattering 
for simplicity. In these approximations, j x is written as 



jx = («A + h)B\ + <J\J\, 
where J x is the mean intensity, defined by 



(7) 



(8) 



Since J x becomes sufficiently close to the Planck function 
in the deep layers, we approximate J x to be B\. Then, 
total energy flux (6) is calculated by 

F\ oc J J J B x [T(x,y,z)]exp(-T Xz )dxdydT Xz 



B\ [T(x,y,z)]exp(-T Xz ) 



x (k x + i x + a x )pdxdydz. 



(9) 

In the following numerical simulations, the integral over 
the disk, dxdy in the above expression, is evaluated as 
summations with respect to the cells, and the integral over 
the beam path, dr Xz , is the summation with respect to the 
thin sliced layers piled along the z-direction. In practice, 
since we are interested in low-degree p-modes with high 
order, the disk image is subdivided into 32 x 32 square 
cells, and the volume is sliced by 200-300 thin zones in 
the z-direction. 

2. 3. Influence of motion 

Our main interest is distortion of the line profiles due to 
oscillations in the stellar atmosphere. Global motion in 
the stellar atmosphere leads to a Dopplcr shift for light 
emitted from there. Besides that, in the case of stellar 
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pulsation, the temperature and any other thermodynami- 
cal quantities vary with the motion, and consequently the 
local surface brightness varies as well. Variations in the 
surface area and in the surface normal due to oscillation 
may also produce some variations of the line profile, but 
these effects other than the Doppler effect are probably 
minor. In this paper, we ignore them for the sake of sim- 
plicity, and restrict ourselves to an investigation of the 
Doppler effect on the line profiles. To take account of the 
Doppler shift corresponding to the line-of-sight velocity, 
we only have to shift the wavelength in equations (4) and 
(9) by the Doppler shift, 



A — > A x [l-v(x,y,z)-e z ]/c, 



(10) 



where v(x,y,z) denotes the velocity field at the coordinate 
(x,y,z), c is the speed of light and e z is the unit vector 
toward the z-direction (see figure 1). The velocity field, 
v, is the sum of the velocity due to oscillation and stellar 
rotation. Convection and turbulence should also be taken 
into account, if any. However, in the case of A-type stars, 
surface convection or turbulence is not expected. 

3. Preparation: Absorption and scattering coef- 
ficients 

3.1. Continuum absorption coefficient 

We restrict ourselves to consider only the wavelength 
range of visual light. In this range, for T > 5000 K, the 
main sources of the continuum absorption coefficient are 
free- free and bound- free absorption of hydrogen atoms and 
of negative hydrogen ions. Following Gray (2005), we cal- 
culate the absorption coefficient in units of area per neu- 
tral hydrogen atom for the free-free absorption and for the 
bound-free absorption, respectively, as follws: 



K(Hff) 



aoA'a^lO"*' 



Ac(H b f ) = ap 



29/ 

A 3 

r> ° 



(11) 



(12) 



where = loge/fcT = 5040/T, / is the ionization energy, 
in units of eV, of hydrogen (13.6 eV), Xn is the excitation 
potential, in units of eV, from the ground state (n = 1) 
to the n-th state [that is, \n = ^(1 — n ~ 2 ))L ancl n o de- 
notes the smallest integer satisfying A < Uq/Roo, where 
i?oo means the Rydberg constant. Here gs and <7bf are 
the Gaunt factors, given by 

0.3456 (XkT 
1+ ,. - „ ,„ - 



gs 



and 



5bf 



1 



(A^oc) 1 / 3 V he 



0.3456 fXR, 



respectively (Menzel 



(13) 



(14) 



k Pckeris 1935), and ao = 
327r 2 e B i? oo /[3 3/2 (47re ) 3 /i 2 c 3 ], where e is the electron 
charge, Sq is the vacuum permittivity, h is Planck's con- 
stant, and k is Boltzmann's constant. 



Table 1. Numerical values of the coefficients {a^} appearing 
in the fitting formula (17). Adopted from Gray (2005). 



coefficient 


numerical value 


a 


1.99654 


«i 


-1.18267 x 10" 5 


n 2 


2.64243 x 10~ 6 


a 3 


-4.40524 x lO" 10 


a>A 


3.23992 x 10~ 14 




-1.39568 x 10~ 14 


a 6 


2.78701 x 10" 23 



As for negative hydrogen, we adopt the following ex- 
pression for the absorption coefficient in units of area (m 2 ) 
per neutral hydrogen atom for the bound-free absorption 
and for the free-free absorption (Gray 2005): 



and 
K(Hr 



4.158 x 10- 31 a bf 



p Q 5/2 10 0.754e 



(15) 



10' 



(16) 



-29p 1 Q/o+/iioge+/ 2 iog^e 

Here, P c is the electron pressure, in units of Pa. The 
wavelength dependence is given for each case in terms of 
a polynomial, which fits sufficiently well to the detailedly 
calculated absorption coefficients (Gray 2005); a b f in the 
case of bound-free absorption given by 



abf 



^2a.iX l , 



(17) 



and fi (i = 0,1,2) in the case of free-free absorption given 
by 



4 

E 

3=0 



fi,3 



log^A, 



(18) 



where the wavelength, A, is given in units of A. The nu- 
merical constants of the coefficients {a,} and the fitting 
formulae for /o, /i, and ji are summarized in table 1 and 
in table 2, respectively. 

It should be noted here that the absorption pro- 
cesses are affected by stimulated emission and the ab- 
sorption actually produced is lower by a factor of {1 — 
cxp[— hc/(XkT)]}. Also, it should be noted here that, fol- 
lowing Gray (2005), the absorption coefficients given in 
this subsection by equations (11), (12), (15) and (16) are 
those for per neutral hydrogen atom. Therefore, we have 
to convert them to those per mass, denoted by k\, dis- 
cussed in the previous section, by multiplying a factor 
nm/p, where nni denotes the number density of neutral 
hydrogen atoms. 

3.2. Line absorption coefficient 

3.2.1. General description 

The mass absorption coefficient, t\, expressed in area per 
unit mass, due to the transition from the Z-th excited level 
of the fc-th ion of the element j, can be written with the 
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Table 2. Coefficients fij of the fitting formulae for /o, /i and fi appearing in the formula (18). Adopted from Gray (2005). 
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1 


2 


3 


4 


fo 


-2.2763 


-1.6850 


0.76661 


-0.053346 




h 


15.2827 


-9.2846 


1.99381 


-0.142631 




h 


-197.789 


190.266 


-67.9775 


10.6913 


-0.625151 



help of the number density of the ions in the state of 
interest, rijki, as 

pix.jki = axnjki {1 - exp [-he/ (A/cT]} 

= a A n J ^^{l-exp[-/ic/(AjfcT)]}, (19) 
rij rijk 

where a\ denotes the absorption coefficient per atom, 
{1 — exp[— hc/(XkT)]} expresses the effect of stimulated 
emission, and rij and rijk are defined by 



and 



E 



n-jki, 



(20) 



(21) 



respectively. Here, njk/nj represents the ionization degree 
of the fc-th ion of the element j, and it is calculated with 
Saha's equation. The ratio njki/njk gives what fraction 
of the ions are at the l-th excited level among the fc-th ion 
of the element j; it is given by 

n jkl _ 9jkl ( Xjkl ' 

kT 



rijk 



U jk (T) 



exp 



(22) 



where gjki denotes the statistical weight of the l-th level, 
Xjkl is the excitation potential of the l-th level, and Ujk(T) 
means the distribution function, defined by 



U jk (T) 



^2g 3 ki exp 



( Xjki\ 
\ kT J 



(23) 



3.2.2. Line broadening 

Each of absorption lines is intrinsically broadened for 
many reasons. One of the main causes is the fact that 
the lifetime of the concerning atomic energy level is not 
infinite due to the uncertainty principle, partly reflecting 
radiation damping, which is often called natural damping, 
due to the interaction of light with atoms, and partly re- 
flecting collisions of atoms absorbing the light with other 
particles. Another main cause is the Doppler effect due 
to the thermal motion of atoms, of which the velocity 
distribution is the Maxwellian distribution. The effects 
of natural damping and of collisional damping lead to a 
Lorentzian profile, characterized by the frequency of the 
line center, z/n, and the damping parameter, T. 

3.2.3. Voigt profile 

The combined effect of the natural damping and the ther- 
mal Doppler broadening is described by a convolution of 
the Lorentzian and the Maxwellian profiles. It is known 
that this convolution integral is reduced to the Hjcrting 
function, H(u,a). which is defined by 



H(u,a) 



t J-oo {u-y) 2 + a 2 



dy. 



(24) 



With the help of this function, the absorption coefficient 
per atom, a\, is expressed as 



1/2 ^of 



-^H{u,a), 



(25) 



47T£om e c 2 Ad 

where Ad denotes the Doppler width corresponding to the 
variance of the Maxwellian velocity distribution of ther- 
mal motion, u = (A — Aq)/Ad, a = TAq/^ttcAd), e and m e 
are electronic charge and mass, respectively, and / de- 
notes the oscillator strength, which is different for each 
atomic level. Most / values are determined from empir- 
ical laboratory measurements. Only in less complicated 
cases, such as hydrogen lines, can theoretical calculations 
be done. 

3.2.4- Damping parameter 

The damping parameter, T, is expressed as the sum of 
that of radiation damping, 7, the quadratic Stark effect, 
74, and the van der Waals effect, 76. As for radiation 
damping, we adopt the classical dipole emission theory 
and set 

8?r e 2 

3 47T£o"lcCAo 

where Ao is the center wavelength of the line. We evaluate 
74 and 76, following Gray (2005), by 



log 74 
and 



18 + ^logC 4 + logP e -^logT 
3 6 



(27) 



log 76 = 19 + 0.41ogC 6 +logP gas - OJlogT 



(28) 

respectively, where P c and P gas are in units of Pa. Here, 
C4 is a constant, but dependent on the line, and Cq is 
given by 

C 6 = 0.3 x 10- 30 [(/ - X ~ Xx)-' 2 -(I- X)- 2 ] , (29) 

where / denotes the ionization energy, \ is the excitation 
potential, and xx = he/ A. /, Xj an d Xx are in units of eV. 
3.2.5. Linear Stark effect for hydrogen lines 

Hydrogen lines are different from other elements' ones 
in the sense that the atomic structure of hydrogen is sub- 
ject to the linear Stark effect. The complete hydrogen 
Stark broadening is the weighted sum of the individual 
Stark components with each showing the Holtsmark dis- 
tribution. We take into account of the linear Stark effect 
in the case of hydrogen lines following Gray (2005). The 
/-values of the unshifted Stark components are denoted as 
/o and the /-values of the shifted components are denoted 
as /±- These values are given in table 3. 
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Table 3. Hydrogen oscillator strengths. Adopted from Gray 
(2005). 



Ha 6562.85 A 



Line 


fo 


f± 


/ 


Ha 


0.248689 


0.392058 


0.640742 


H/3 


0.00 


0.119321 


0.119321 


H 7 


0.007014 


0.037656 


0.044670 


H<5 


0.00 


0.022093 


0.022093 



The approach of Stehle (1994) is used for H/3, while 
that of Vidal ct al. (1973) is used for the other Balmer 
lines. With the help of these Stark profiles, we calculate 
the absorption coefficient for hydrogen by 



"A 



p 2 A 2 
7T 

4:7Teom c c 2 Ad 



S(AX/E ) 



H(u,a) + f H(u,a) 



where * means convolution. Here, 

2/3 



E, 



( Aim 



(30) 



(31) 



\ 3 J 47T£o ' 

where n denotes the number density of ions. 
3.3. Scattering coefficient 

The contribution of Rayleigh scattering is substantially 
large at shorter wavelengths ,such as A < 300 nm. But, 
for the optical region, it is negligibly small. We therefore 
only have to consider Thomson scattering by free elec- 
trons. The scattering coefficient is then independent of 
wavelength, and is given by 

a A =^ ( g2 V ne , (32) 
3 \47re w e c 2 / p 

where n c denotes the electron number density. 

4. Check of our numerical code 

4-1. Intrinsic line profiles 

In order to check our numerical code, we calculated first 
the line profile of the Balmer lines by using a solar model 
(Christcnscn-Dalsgaard et al. 1996). The upper panel 
of figure 2 shows the computed line profile of H a line. 
In this procedure, the effect of rotational broadening was 
not taken into account, and micro-turbulence was ignored. 
The values of 74 and 76 were adjusted so that the wing 
part of the computed profile fit well to the observed one 
(e.g., Babcock k Moore (1947)). We confirmed how 
well the computational programs developed here do work, 
though there remains some discrepancy. Fitting is sub- 
stantially improved in the case of H/3 and H7. 

How much does each layer of the atmosphere contribute 
to formation of the spectral line? This can be seen by 
the integrand of the integral with respect to t\ z , near 
the disk center, in equation (9). It is often called the 
contribution function. The lower panel of figure 2 shows 
the contribution functions for H a line for the solar model. 
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Fig. 2. Upper: Computed line profile of the Ha line in a 
wide range. The ordinate is the flux normalized with the 
continuum flux, and the abscissa is the wavelength measured 
from the line center, in units of A. Lower: Contribution func- 
tions for H a line for the solar model. These are the integrand 
of the integral with respect to t\ z in equation (9), evaluated 
at the wavelengths, which are marked along the computed 
line profile by dots in the upper panel. The abscissa is the 
optical depth for A = 5000 A, rsooo- 

These are the integrand of the integral with respect to t\ z 
in equation (9), evaluated at the wavelengths, which are 
marked along the computed line profile by dots in the 
upper panel of figure 2. The abscissa is the optical depth 
for A = 5000 A, 

r 5000- The contribution function for the 
line center concentrates around T5000 — 10~ 3 ~~ 4 , which 
corresponds to 300 ~ 500 km above the photosphere. 

4-2. The LPV caused by radial pulsation 

Let us see the influence of pulsation on the line profile of 
Ha line thus calculated. The current purpose is to see 
how much difference in the line profile variation is caused 
by the difference in the vertical wavelength of motion. For 
this reason, we assume radial pulsation of a non-rotating 
star, in which case no lateral motion exists. 

We consider three cases: (i) a standing wave with a con- 
stant pulsation amplitude, which is regarded as represent- 
ing a case of very long vertical wavelength; (ii) a standing 
wave with the wavelength of the order of the thickness 
of the line forming layer; (iii) an upwardly running wave 



No. 
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Standing wave 
(Long wavelength) 



Standing wave 
(Short wavelength) 



Runing wave 
(Short wavelength) 



1.0 




AX [A] AX [A] AX [A] 

Fig. 3. Theoretically expected line-profile variation (LPV) of Ha line produced by a sinusoidally oscillating radial mode of a 
nonrotating star. Left: Pulsation velocity amplitude was assumed to be spatially constant. This is regarded as representing a case 
of very long vertical wavelength. Middle: Oscillation was assumed to be a standing wave with the wavelength of the order of the 
thickness of the line forming layer. Right: Oscillation was assumed to be an upwardly running wave with the wavelength of the order 
on the thickness of the line forming layer. In each panel, the differences between the individual profiles at a sequence of 100 phases 
and their mean are stacked with phase increasing downwards. A short bar in the left panel indicates the length corresponding to 
10% of the continuum intensity. Above them are their mean (dashed line) and the line profile at the phase zero (solid line), and 
below them is plotted their standard deviation. The abscissa is the wavelength from the line center of the intrinsic profile, in units 
of A. The left-hand ordinate scales refer to the phase of the oscillation. 



with the wavelength of the order of the thickness of the 
line forming layer. 

The line profile variation is most clearly seen in the 
time series of the residuals that remain after subtracting 
from each line profile its temporal mean. In the left panel 
of figure 3 is displayed the line-profile variation (LPV) of 
the case (i) , obtained from the current calculation using a 
solar model. Plotted is the deviation of the integrated ra- 



diant intensity from its mean as a function of wavelength 
for a sequence of values of the pulsation phase extend- 
ing over one complete period. Also, the line profile at 
the phase zero is displayed to compare with the mean 
profile averaged over one period. Since the pulsation am- 
plitude is assumed to be spatially constant in this case, 
the line-forming layer moves up and down as a whole. 
Consequently, the spectral line shifts from blue to red and 
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Standing wave Standing wave Runing wave 

( Long wavelength ) ( Short wavelength ) ( Short wavelength ) 
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Fig. 4. Theoretically expected variation in the bisector of line profiles during a one complete period. In each panel, the line profiles 
(solid lines) and the bisectors (dashed lines) at the pulsation phases 0.0, 0.1, 1.0 are stacked with phase increasing downwards. 
The abscissa is the deviation from the line center of the intrinsic profile, and the ordinate of each diagram at a fixed pulsation phase 
is the line intensity normalized with the continuum level. Left: Case of a constant velocity amplitude. Middle: Case of a standing 
wave with the wavelength of the order of the thickness of the line forming layer. Right: Case of an upwardly running wave with the 
wavelength of the order of the thickness of the line forming layer. A short bar in the left top panel indicates the length corresponding 
to 20 % of the continuum intensity. 



back again, while keeping its symmetric profile. This be- 
havior can be clearly seen in the motion of the bisector 
of the line profile, which is displayed in the left panel of 
figure 4. The bisector is almost always a straight line, and 



it simply moves from side to side with time . 

The small deviation from a straight line reflects the spherical 
effect of stratification. The angle between the normal to an at- 
mospheric layer and the line of sight gradually becomes slightly 
larger with height of the layer. Hence, the Doppler shift is 
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Fig. 5. Variation in the equivalent width of Ha lino in the 
case that lateral distribution of hydrogen in the atmosphere 
is proportional to |cos<9|, where 9 denotes the colatitude with 
respect to the magnetic axis, which is assumed to be perpen- 
dicular to the rotation axis assumed to be perpendicular as 
well to the line of sight (i = 90° and = 90°). The rotation 
phase zero is the phase when the magnetic axis is oriented 
toward the line of sight. The ordinate is in units of A. 

The middle panel of figure 3 displays the LPV seen in 
the case of (ii) . The difference from the case (i) would be 
hardly perceptible in this form, but it is clearly seen in 
the motion of the bisector, which is shown in the middle 
panel of figure 4. The vertical wavelength of pulsation has 
been assumed in this case to be as short as the thickness 
of the line-forming layer. As a result, the velocity varies 
with height, and the line shape no longer has symmetry. 
Furthermore, there is a nodal surface at a certain level. 
These features are well recognized in the middle panel of 
figure 4. The right panel of figure 3 displays the LPV seen 
in the case of a running wave, propagating upward in the 
atmosphere. The running wave feature is clearly seen in 
the motion of bisector (the right panel of figure 4) . On the 
contrary to the case of a standing wave, the line depth at 
which the bisector crosses the line center of the intrinsic 
line profile docs not stay still, but shifts downward with 
time. 

4-3. Influence of lateral chemical inhomogeneity and ro- 
tation of the star 

If the chemical distribution in the stellar atmosphere is 
laterally homogeneous, the effect of stellar rotation on a 
stellar spectrum is to broaden spectral lines. The amount 
of broadening depends on the rotation velocity and the the 
inclination of the rotation axis to the line of sight, while 
the equivalent width is kept independent of rotation. As 
explained in section 1, the atmospheres of magnetic Ap 
stars are laterally inhomogeneous and some chemical ele- 
ments are selectively levitated near the magnetic polar re- 
gion. To see the effect of such lateral inhomogeneity upon 
LPV, we assume here that hydrogen is slightly concen- 
trated in the magnetic polar region. The number density 

slightly larger with line intensity even if the pulsation velocity 
amplitude is a constant. 



of hydrogen is supposed to vary proportionally to |cos0|, 
where 8 is the colatitude with respect to the magnetic 
axis. Naturally enough, the equivalent width of H a line 
is dependent on the aspect angle. As a consequence, if the 
magnetic axis is misaligned with the rotation axis of the 
star, the equivalent width of hydrogen lines synchronously 
varies with the rotation. To see the largest effect, we 
demonstrate here, in figure 5, the most extreme case, in 
which the line opacity is proportional to |cos#|, though it 
is implausible in the case of H a line. 




y , ▼ 

-8-6-4-2 
A?. [A] 



1.0 




-2.0 -1.0 0.0 1.0 2.0 

AA. [A] 



Fig. 6. Theoretically expected LPV of Ha line in the case 
that hydrogen is slightly concentrated near the magnetic polar 
region in the atmosphere. Upper: The individual profiles at 
a sequence of 100 phases with phase increasing downwards. 
Lower: The differences between the individual profiles at a 
sequence of 100 phases, and their mean are stacked with phase 
increasing downwards; above them are their mean (dashed 
line) and the line profile at the phase zero (solid line), and 
below them is plotted their standard deviation. In each panel, 
the abscissa is in units of A, and the left-hand ordinate scales 
refer to the rotation phase. A short bar in each panel indicates 
the length corresponding to 1 % of the continuum intensity. 
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Figure 6 shows the line profiles in the case of i = 90° 
and (3 = 90° at a sequence of 100 rotation phases (the 
upper panel) and the differences between the individual 
profiles and their mean (the lower panel). When one 
magnetic polar region is coming to appear on the visi- 
ble hemisphere while the other polar region is coming to 
disappear from there, the apparent distribution of the el- 
ement over the visible hemisphere becomes conspicuously 
inhomogeneous. As a consequence, the line shape becomes 
asymmetric at such a rotational phase. 

5. Calculation of line-profile variation in A-type 
stars 

5.1. Equilibrium model for A-type stars 

Based on the numerical technique explained in the pre- 
vious sections, let us proceed to the case of Ap stars. 
We choose a 2 M star at the zero age as an equilibrium 
model. The star is assumed to be chemically homoge- 
neous throughout the star, and the hydrogen mass ratio, 
X, and the metal mass ratio, Z, are assumed to be 0.73 
and 0.02, respectively, and the chemical composition of 
the solar abundance is assumed. The equilibrium model 
of the star is constructed with a program originally written 
by Paczynski (1970). The opacity table in the program 
is, however, updated, and we adopt the OPAL (Iglesias et 
al. 1996) for a wide range of temperature and supplement 
the low temperature range (logT < 4.0) with Ferguson et 
al. (2005). The global physical parameters of the model 
are summarized in table 4. 

Table 4. Physical properties of the equilibrium model. 

M/Mq R/Rq L/L q T cff X Z 
2.0 1.64 14.7 8800K 0.73 0.02 



5.2. Ionization degrees of hydrogen and neodymium ions 

It should be noted here that, contrary to the solar case, 
the density in the high atmosphere is so low that hydro- 
gen is ionized there (r 50 oo > 10~ 5 ), as shown in figure 8. 
Consequently, the H a line in A-type stars is more sensi- 
tive to the very outer layer than in the case of the Sun. 
This is clearly seen in the contribution function shown in 
figure 9, which should be compared with figure 2, display- 
ing the contribution function for the solar model. 



3000 




T 5000 



Although the atmospheric structure of Ap stars has not 
been definitely established, it may be substantially deviate 
from the standard one due to the strong magnetic fields 
and peculiar chemical compositions. Following Shibahashi 
k Saio (1985), we have adopted an analytic form of the T- 
t relation for the atmosphere, which was fitted to Kurucz's 
(1979) model atmosphere: 

T 4 = ^Tcff [r + 0.9052 - 0.3367exp(-2.54r) 

-0.1645cxp(-300r) - 0.1440cxp(-10000r). (33) 

The top panel of figure 7 displays the relation between the 
optical depth at A = 5000 A and the geometrical depth 
from the top layer for the present model. 

The sound speed profile in the atmosphere is displayed 
in the bottom panel of figure 7. The sound speed is about 
lOkms -1 at the photosphere, and becomes as slow as 
7kms~ 1 in the high atmosphere. In the case of acous- 
tic modes with a period of ~ 600 s, the wavelength is as 
short as the thickness of the atmosphere. 



Fig. 7. Top: Depth from the top layer of the present model 
as a function of the optical depth at A = 5000 A. Middle: 
Density profile in the atmosphere as a function of the optical 
depth at A = 5000 A. Bottom: Sound speed in the atmosphere 
as a function of the optical depth at A = 5000 A. 

Recent high-speed spectroscopic observations of oscil- 
lations of roAp stars have revealed that spectacular line- 
profile variation are seen in the spectrum lines of Pr III 
and Ndlll. It is these observations that have motivated 
the work reported in this paper. We then pay atten- 
tion to Ndlll. Figure 8 displays the ionization degrees 
of neodymium ions, as well as that of hydrogen, in the 
present A-star model. The ionization degree of Nd III goes 
up with height after reaching a local minimum at a few 
hundred kilometers above the photosphere, and almost all 
of neodymium is Nd III in the high atmosphere. This is 
the reason why Ndlll lines show much spectacular line- 
profile variation than Ndll lines. In these calculations, 
the atomic data of Nd were taken from a database offered 
by US National Institute of Standards and Technology 
(http://physics.nist.gov). However, the statistical weight 
of Nd IV was not available there. We took the value from 
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and some stars show overabundance of 4 dex (Adclman 
1973). In many cases, the observed Ndlll lines with low 
excitation energies are strong, though Nd II lines are much 
weaker. This implies that neodymium is mainly concen- 
trated in a thin layer above T5000 < 1CP 5 . Ryabchikova 
et al. (2002) modeled the spectra of one of Ap stars, 
7Equ, by supposing log(ANd/An) — —4 for logr < —9 
under the assumption of LTE, where iV^d/^Vtotai denotes 
the fraction of neodymium. Mashonkina, Ryabchikova 
& Ryabtsev (2005) modeled the same star with the as- 
sumption of NLTE by supposing log^Nd/A^ii) — ~~ 7 for 
logr < -4. 



Ndlll 6145.07 A 



Fig. 8. Ionization degrees of neodymium and of hydrogen in 
the present model of an A-star. Plotted as a function of the 
optical depth tsooo- 
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Fig. 9. Upper: Computed line profile of the Hct line in a 
wide range. The abscissa is the wavelength measured from 
the line center, in units of A. Lower: Contribution functions 
for Ho line for present model of an A-star. These are the 
integrand of the integral with respect to t\ z in equation (9), 
evaluated at the wavelengths, which are marked along the 
computed line profile by dots in the upper panel. The abscissa 
is the optical depth for A = 5000 A, Tsooo- 



Kurucz's ATLAS and set logg ; 
5.3. Contribution functions 



1.8 for any energy levels. 



The abundance of neodymium involved in Ap stars seems 
to generally be much larger than the solar case, in which 
log(AT Nd /Ar H ) + 12 = 1.42 ± 0.04 (Asplund ct al. 2009), 
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Fig. 10. Upper: Line profile of Ndlll A 6145.07 A line for 
the current model with log(A r Nd/^ r H) = — 8- The abscissa is 
the wavelength measured from the line center, in units of A. 
Lower: Contribution functions for the above line. These are 
the integrand of the integral with respect to t\ z in equation 
(9) , evaluated at the wavelengths marked along the computed 
line profile by dots in the upper panel. The abscissa is the 
optical depth for A = 5000 A, tsooo- 

In this paper, we assume, however, the vertically uni- 
form distribution of elements for the sake of simplicity, 
for our current main purpose is to demonstrate the use- 
fulness of analyses of LPV as an information source of os- 
cillation and atmospheric structures of roAp stars, rather 
than to model the individual star in detail. The contri- 
bution functions for Ndlll A6145.07A line in the case of 
log^Nd/A^H) = —8 and in the case of log^Nd/A^i) = — 5 
are shown in figures 10 and 11, respectively. The oscilla- 
tor strength was taken from Cowley & Bord (1998). The 
effect of rotation of the star was ignored in these calcula- 
tions. 
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Fig. 11. Same as figure 10, 



but for the case of 



5.4- Line-profile variation due to nonradial oscillations 
in A-type stars 

As explained in section 1, the general features of pulsa- 
tions in roAp stars are essentially well explained in terms 
of high order, axisymmetric dipolc p-modes, of which the 
symmetry axis is aligned to the magnetic axis of the star, 
which is inclined to the rotation axis of the star. In this 
subsection, we simulate line-profile variation of the H a 
line and Ndlll A 6145.07 A due to such oscillation modes 
in roAp stars by using the equilibrium model described in 
the previous subsection. Both of the elements are assumed 
to slightly concentrate in the magnetic polar region, dis- 
tributed varying proportionally to cos 2 #, where 9 denotes 
the polar angle measured from the magnetic axis. The 
following simulations are for the case of line opacity being 
proportional to cos 2 (9 + 1/2. We treat the inclination angle 
between the rotation axis of the star and the symmetry 
axis of pulsation as a free parameter, /?. Also the angle 
between the line-of-sight and the rotation axis of the star 
is treated as another free parameter, i. The aspect angle 
of the symmetry axis of pulsation, a, gradually changes 
with time due to rotation of the star, as follows: 



cosa = cosicos/3 + sin i sin /3 cos fit, 



(34) 



where f2 denotes the rotational frequency. We will show 
the theoretically calculated line profiles at every 0.125 ro- 
tation phase. 



5.^.1. Standing eigenmode 

We consider here a standing eigenoscillation of the A-star 
model. For the sake of simplicity, the oscillation is as- 
sumed to be adiabatic. Also, we ignore the effects of ro- 
tation of the star and treat the star as being spherically 
symmetric, and deal with only low-degree, high-order p- 
modes, of which period is 600 s; — in a typical period range 
of roAp stars. In this case, the Eulerian perturbation to 
the gravitational potential is negligibly small, and then 
the basic equations governing the linear, adaiabatic oscil- 
lations of a star are well described with only equations 
of conservations of momentum, mass, and entropy; they 
are reduced to a form expressed with only the velocity 
perturbation, v; 



d 2 v 



V(v ■ Vp) ■ 
1 



V(/oc 2 V •«) 



- - (v ■ V p)Vp - (V ■ vWp, 
P 



(35) 



where p and p denote the pressure and the density, respec- 
tively, and c means the sound speed. We solve numerically 
equation (35) with a set of proper boundary conditions 
to compute the eigenfunctions of linear, adiabatic, dipole 
oscillations of the A-star model, following the procedure 
described in Unno et al. (1989). The amplitude is arbi- 
trarily chosen in the linear calculation. We choose it so 
that the maximum velocity is still subsonic; \v \ = 7kms~ 1 
at the top layer of the atmosphere. 

In the outer envelope of a star, the effect of sphericity 
becomes unimportant. To understand the oscillation fea- 
ture there, it is instructive to consider a simple solution 
of plane waves of the form 



v oc exp 



2H r 



exp(ik ■ x + iu)t) 



(36) 



with an assumption of k z 3> H p , where H p = —dz/dhip is 
the density scale height and lu is the frequency. The expo- 
nentially growing factor with height in equation (36) arises 
so as to conserve wave energy in vertical direction since the 
density in atmosphere decreases with height. Substituting 
equation (36) into equation (35), we obtain a local disper- 
sion relation, which gives the vertical wavenumber k z for a 
given set of the horizontal wavenumber and the frequency. 
For high-order p-modes, the dispersion relation is reduced 
to 



where 



\z) {c 



?(*)}: 



4# 2 



1 



2 dHp 

dz 



(37) 



(38) 



denotes the squared critical frequency. Equation (37) 
means that the oscillation is evanescent and the wave en- 
ergy, p\v\ 2 , decreases monotonically with height if u 2 < 
w 2 (z). The latter condition is required for the wave to 
be standing wave. Since what we are now considering are 
standing cigenmodes, they are trapped waves that are re- 
flected at the two boundaries of the system. Hence, they 
become evanescent above the photosphere, and then have 
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Fig. 12. Theoretically expected LPV of the Ha line produced by an axisymmetric dipole (Z = l,m = 0) mode at a sequence of the 
rotation phase </>=0.0, 0.125, 0.25, 0.375, and 0.5 (from left to right). The rotation axis is assumed to be inclined with respect to the 
linc-of-sight by 90° , and the symmetry axis of pulsation is assumed to be inclined with respect to the rotation axis of the star by 
90°. In each panel, the differences between the individual profiles at a sequence of 100 pulsation phases and their mean are stacked 
with phase increasing downwards; above them is their mean, and below them is plotted their standard deviation. The abscissa is 
the wavelength measured from the line center in units of A. The left-hand ordinate scales refer to the phase of the oscillation, in 
units of the pulsation period. A short bar in the left panel indicates the length corresponding to 1 % of the continuum intensity. 



no node in the line-forming layer. However, the frequen- 
cies are so close to the critical frequency near the surface 
that the e-folding scale of attenuation of evanescent wave 
energy is quite long. Hence, the wave energy decreases 
only very gradually with height. As a consequence, the 
velocity amplitude increases with height because the de- 
crease in density is steeper than the energy density (see 
middle panel of figure 7). 

It should be noted that the presence of a node in the 
line-forming layer has been observationally shown in some 
roAp stars (see Sect.l and 6.2). Also it should be noted 
that, by taking account of magnetic effects on pulsation 
and an atmospheric model having a temperature inversion 
layer, Saio, Ryabchikova & Sachkov (2010) demonstrated 
magneto-acoustic modes might have a node layer in the 
line-forming layer. 

Figure 12 displays the line-profile variation of the Ha 
line produced by an axisymmetric dipole p-modc at a se- 
quence of the rotation phase 4> =0.0, 0.125, 0.25, 0.375 
and 0.5 (from left to right). The rotational velocity at the 
equator is assumed to be 8.3 krns" 1 . The rotation axis 
has been assumed to be inclined with respect to the line- 
of-sight by 90° , and the symmetry axis of pulsation has 
been assumed to be inclined with respect to the rotation 
axis of the star by 90°. 

Figure 13 displays the line-profile variation of Ndlll 
A6145.07 A line produced by the same mode and the same 



geometrical configuration as in the case of figure 12 at a 
sequence of the rotation phase <fi =0.0, 0.125, 0.25, 0.375 
and 0.5 (from left to right). As seen in this figure, the 
residual is not symmetric at some rotational phases. This 
is caused by lateral inhomogeneity of the distribution of 
neodymium. Since we have assumed that the element is 
slightly concentrated in the magnetic polar region, which 
is misaligned with the rotation axis of the star, the line 
shape becomes significantly asymmetric when one mag- 
netic polar region is coming to appear on the visible hemi- 
sphere while the other polar region is coming to disap- 
pear. The same effect can be seen in the Ha line shown 
in figure 12. This effect is most conspicuous in the case of 
i = 90° and /3 = 90°. 
5.4-2. Running wave 

Though the cigenmodes of the present model arc evanes- 
cent above the photosphere and have no node there, spec- 
troscopic observations of rapid oscillations of roAp stars 
indicate the presence of a nodal surface in the line-forming 
layers (e.g., Baldry et al. 1998, 1999). Furthermore, in 
some stars, while Ha line shows the standing wave fea- 
ture, the absorption lines of Ndlll and some other rare 
earth elements imply the running wave feature (Kurtz, 
Elkin & Mathys 2003). Our present model is not fit to 
these Ap stars. The presence of a strong magnetic field 
and/or peculiar chemical composition in the atmosphere 
is likely to change substantially the atmospheric struc- 



14 



= 0.000 




T. Nomura, J. Naito & H. Shibahashi [Vol. , 

= 0.125 = 0.250 = 0.375 = 0.500 




0.0 
A*. [A] 




Fig. 13. Same as figure 12, but for the Ndm A6145.07A line. 



ture, and particularly the oscillation feature. In the case 
of a star with B = 10 mT (= lkG), the local Alfven speed 
above the photosphere exceeds the sound speed. As a 
consequence, the acoustic wave and the magnetic wave 
couple with each other around the photosphcric level, and 
the eigenmode in such a system is expected to be a mode 
with energy leakage (Sousa & Cunha 2008; Shibahashi 
2008). We consider here an outwardly running wave and 
see its effect on the spectral line profile. 

We suppose that the wave energy is conserved during 
propagation through the line-forming layer. Then, the 
velocity amplitude grows with height in inverse proportion 
to the square root of density. Though, in the very outer 
layer, the amplitude may become supersonic, we assume 
here that the amplitude is still subsonic and sinusoidal, 
and estimate the vertical wavelength from the frequency 
and the sound speed. Figure 14 displays the line-profile 
variation of the H a line produced by such a running wave 
at a sequence of the rotation phase 4> =0.0, 0.125, 0.25, 
0.375 and 0.5. The angle between the rotation axis of 
the star and the line-of-sight has been assumed to be 90°, 
and the symmetry axis of pulsation has been assumed 
to be inclined with respect to the rotation axis of the 
star by 90°. With propagation of the wave, the layer of 
the maximum amplitude shifts upward. Comparied with 
the line-profile variation in the case of a standing wave 
shown in figure 12, we see that the running wave feature 
is characterized by the ridge structure migrating from blue 
to red and then back again to blue. 

The Ndm A 6145.07 A line forms at higher levels than 
the H a line, and so its line profile is deformed conspicu- 
ously around the phase at which the amplitude becomes 



large at the high atmosphere. This is seen as a spotty pat- 
tern shown in figure 15, in which case we assume i = 90° 
and /3 = 90°. 

6. Variation in bisectors 

6.1. Simulation results 

The variation in the bisector on the pulsation phase has 
been measured in some roAp stars. Let us consider the bi- 
sector motion expected from the present simulation. The 
most left panel of figure 16 shows the variation in the bi- 
sector of the H a line caused by the standing wave of ax- 
isymmetric, dipole, eigenmode on a sequence of pulsation 
phases in the case of i = 90° and (3 = 90°. The ampli- 
tude of variation of the bisector gradually becomes larger 
with decrease in the residual intensity. The largest ampli- 
tude of variation in the bisector is apparently substantially 
smaller than the assumed largest velocity amplitude in the 
top of the atmosphere. This is partly because the density 
in the line-forming layer is higher than at the top of the 
atmosphere, and partly because the Dopplcr effect reflects 
only the line-of-sight velocity. As seen in the second left 
panel of figure 16, the variation in the bisector in the case 
of a running wave is even smaller than in the case of a 
standing wave. This is because the standard deviation of 
the amplitude at a fixed height is smaller than that in the 
case of a standing wave. 

The variation in bisector of Nd III A6145.07 A line on the 
pulsation phase is shown in figure 16. The second-right 
panel is the case of the standing wave of axisymmetric, 
dipole mode. The geometrical configuration is the same 
as in the case of the most-left panel. As seen in this panel, 
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Fig. 14. Theoretically expected LPV of the H a line produced by an outward running wave of the axisymmetric dipole (Z = l,m = 0) 
oscillation at a sequence of the rotation phase <j> = 0.0, 0.125, 0.25, 0.375, and 0.5 (from left to right). The rotation axis is assumed 
to be inclined with respect to the line-of-sight by 90°, and the symmetry axis of pulsation is assumed to be inclined with respect 
to the rotation axis of the star by 90°. In each panel, the differences between the individual profiles at a sequence of 100 pulsation 
phases and their mean are stacked with phase increasing downwards; above them is their mean, and below them is plotted their 
standard deviation. The abscissa is the wavelength measured from the line center in units of A. The left-hand ordinate scales refer 
to the phase of the oscillation, in units of the pulsation period. A short bar in the left panel indicates the length corresponding to 
0.1% of the continuum intensity. 



the general tendency is similar to the case of the H a line, 
though the range of variation is slightly smaller. 

The variation in the case of a running wave is displayed 
in the most-right panel. Compared with the case of the 
Ha line, the amplitude of variation becomes more con- 
spicuous with a decrease in the residual intensity. This is 
because the line center of Ndlll A6145.07A line is formed 
around T5000 — 10 -65 , while the line center of Ha line is 
sensitive to the layer around r 50 oo — 10~ 5 ' 5 . 

6.2. Comparison with observations 

Mkrtichian, Hatzes & Kanaan (2003) found in the roAp 
star HD 137949 (33 Lib) that the oscillation phases of 
Ndlll lines arc nearly 7r out of phase with Ndll lines. 
This indicates presence of a nodal layer just beneath the 
Ndlll line-forming layer of the star. They also suspected 
the presence of another nodal layer close to the contin- 
uum formation level. Similar results were obtained by 
Kurtz, Elkin & Mathys (2005a). They found that there 
is a phase jump between the Nd II line- forming layer and 
the Nd III line-forming layer, and concluded the presence 
of a nodal layer between these two types of lines. They 
also found that the amplitude of the Ndlll line, which is 
formed in the higher layer than the Ndll line, is almost 
constant, although a small increase with height cannot 
be ruled out. The present calculation, however, does not 



show any node just beneath the Ndlll line-forming layer, 
and this obviously means that fine tuning in making the 
equilibrium models and in mode calculation is necessary. 
Naively, an increase of the amplitude with height from 
the node is theoretically expected. Furthermore, Kurtz, 
Elkin & Mathys (2005a) showed that the phase of the 
Ndlll line is constant in the high atmosphere, but shows 
an outwardly running wave character in the deeper layer. 
This implies that the magnetic effect on the oscillations 
should be taken into account. 

Elkin, Kurtz & Mathys (2005) found in the case of the 
roAp star HD 99563 that the amplitude and phase of H a 
core vary strongly with the line depth. There seems a 
nodal layer within the Ha line-forming layer. The same 
feature was found in some other rare erath lines including 
Ndlll A6145.07A line. The present model does not show 
a node in the Ha line-forming layer. This again means 
that fine tuning in making the equilibrium models and in 
mode calculation is necessary. 

Ryabchikova et al. (2007a) analyzed the bisector mo- 
tion of the H a and neodymium lines of some roAp stars. 
Qualitatively, the general feature of the bisector motion of 
H a line is in good agreement between the simulation and 
the observations, in the sense that the amplitude grad- 
ually becomes larger with the line depth. However, the 
situation is quite different in the case of the Ndlll line. 
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Fig. 15. Same as figure 14, but for the Ndlll A6145.07 A line and the angle i = 90° and f3 = 90° 



While the present simulation shows that the amplitude 
becomes larger with the line depth, the observations show 
the opposite tendency. This implies a possibility of strong 
dissipation in the high atmosphere, but more careful com- 
parison should be made before concluding. 

Kurtz et al. (2007) found an increase in amplitude of 
the Nd III lines with line strength in the case of HD 134214, 
and also found phase variations of those lines, suggesting 
an outwardly running wave. 

Elkin, Kurtz & Mathys (2008) studied LPV of the roAp 
star HD 176232 (10 Aql), and showed that the amplitude 
of the H a line increases with the atmospheric height, and 
that the phase variation on the atmospheric height indi- 
cates an outwardly running wave feature. They also found 
a similar tendency concerning the amplitude of Nd III lines 
to that found by Ryabchikova ct al. (2007a). That is, 
they found that the phase of the weak lines is almost 
constant, suggesting a standing wave, and that the am- 
plitude decreases with the atmospheric height. Lines of 
Nd III with different intensities were found to have pulsa- 
tion amplitudes different from the weaker lines that are 
formed more deeply in the atmosphere having higher am- 
plitudes. It was also found that the bisectors for strong 
Ndlll lines show significant changes of phase, indicating 
complex variations in the oscillation phase as a function 
of the atmospheric height. 

7. Discussion 

A simulation of LPV beyond single surface treatment is 
definitely necessary to obtain quantitative information of 
roAp stars from spectroscopic observations quoted in the 



previous section 6.2 with high time resolution, high spec- 
tral dispersion, and a high signal-to-noise ratio. This is 
a new information source of Ap stars, and is expected to 
be a promising new tool of asteroseismology to probe 3-D 
structure of atmospheres of these enigmatic stars. Our 
aim is to establish an inversion approach of LPV to probe 
Ap stars' atmosphere and pulsation. This present paper 
is the very first step for this purpose, and there is much 
room for improvement. The following is a road map to 
our goal. 

• Levitation of chemical elements due to radiative 
pressure should be taken into account in model at- 
mosphere. 

• The Zecman effect should be taken into account in 
calculating the line shapes of spectral lines. 

• Nonadiabatic mode calculation should be done in- 
stead of adiabatic calculations. 

• Magnetic effects on oscillations should be taken into 
account. 

• The effect of temperature perturbation should be 
taken into account in LPV calculations. 

• More systematic calculations should be carried out 
for many more lines. 
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Fig. 16. Theoretically expected variation in bisector of line profiles caused by axisymmteric dipole oscillation during a one complete 
period. The rotation axis and the pulsation axis are assumed to be aligned with the line-of-sight (i = /3 = 90°). From left to right: Ha 
in the case of a standing eigenmode; Ho in the case of a running wave; Ndm A 6145.07 A line in the case of a standing eigenmode; 
Ndm A 6145.07 A line in the case of a running wave. In each panel, the line profiles (solid lines) and the bisectors (dashed lines) at 
the pulsation phases 0.0, 0.1, 1.0 are stacked with phase increasing downwards. The abscissa is the deviation from the line center 
of the intrinsic profile, in units of A, and the ordinate of each diagram at a fixed pulsation phase is the line intensity normalized 
with the continuum level. A short bar in the left top panel indicates the length corresponding to 20% of the continuum intensity. 
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